##############################################################################
## bacom
##
## Description:
##   'bacom' calls Java class Noticcor to perform normal tissue contamination
##   correction using Bayesian analysis of copy number mixtures (BACOM).
##   'bacom' requires R package 'rJava' and Java Virtual Machine installed
##   and configurated properly.
## 
## Usage:
##   alpha <- bacom(chrID, copyNumber, alleleA, alleleB, isGenotypeAB, segMean)
##
## Arguments:
##   chrID:        a character vector of chromosome ID, the same length as
##                 copyNumber
##   copyNumber:   a numeric vector of DNA copy numbers (the normal segments
##                 have mean 2)
##   alleleA:      a numeric vector of allele A, the same length as copyNumber
##   alleleB:      a numeric vector of allele B, the same length as copyNumber
##   isGenotypeAB: a boolean vector indicating whether the corresponding loci
#                  are genotype AB, the same length as copyNumber
##   segMean:      a numeric vector of segment means, the same length as
#                  copyNumber
##
## Value:
##   'bacom' returns the estimated fraction of normal cells in the sample.
##
## Example:
##   alpha <- bacom(chrID, copyNumber, alleleA, alleleB, isGenotypeAB, segMean)
##
## Author:
##   Bai Zhang    baizhang@vt.edu
##
##############################################################################


bacom <- function(chrID, copyNumber, alleleA, alleleB, isGenotypeAB, segMean) {

  if (!library("rJava", logical.return=TRUE)) {

    cat("\n", "  Failed to load package 'rJava'. \n\n")
    return(NULL)
    
  }

  if (.jinit(classpath = "." ) < 0) {

    cat("\n", "  Failed to initialize the Java Virtual Machine. \n\n")
    return(NULL)

  }

  chrID <- paste(chrID)   ## convert to character vector, if it is numeric
  
  jbacom <- .jnew("bacom/Noticcor", .jarray(chrID), .jarray(copyNumber),
                  .jarray(alleleA), .jarray(alleleB),
                  .jarray(isGenotypeAB), .jarray(segMean))

  if (!is.jnull(jbacom)) {

    .jcall(jbacom,, "correctNTC")
    alpha <- .jcall(jbacom, "D", "getNormalTissueFraction")

  } else {
    cat("\n", "  Failed to create Java object of class Noticcor. \n\n")
    return(NULL)
  }


  segment.alpha <- numeric(0)
  del.seg.index <- 1

  while (jbacom$getDeletionSegNormalTissueFraction(as.integer(del.seg.index)) > 0 || is.nan(jbacom$getDeletionSegNormalTissueFraction(as.integer(del.seg.index)))) {
    segment.alpha <- c(segment.alpha, jbacom$getDeletionSegNormalTissueFraction(as.integer(del.seg.index)))
    del.seg.index <- del.seg.index + 1
  }


  for (i in 1:length(segment.alpha)) {
    cat("\nDeletion segment ", i, " alpha = ", round(segment.alpha[i], 3), "\n")
  }
  
  cat("\nalpha = ", alpha, "\n\n")
  
  alpha
  
}


